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Abstract 

(**^ The optimal control of a mechanical system is of crucial importance in many realms. Typical examples 

are the determination of a time-minimal path in vehicle dynamics, a minimal energy trajectory in space 

amission design, or optimal motion sequences in robotics and biomechanics. In most cases, some sort 
of discretization of the original, infinite-dimensional optimization problem has to be performed in order 
to make the problem amenable to computations. The approach proposed in this paper is to directly 
I discretize the variational description of the system's motion. The resulting optimization algorithm lets 

the discrete solution directly inherit characteristic structural properties from the continuous one like 
symmetries and integrals of the motion. We show that the DMOC approach is equivalent to a finite 
OO difference discretization of Hamilton's equations by a symplectic partitioned Runge-Kutta scheme and 

C**) employ this fact in order to give a proof of convergence. The numerical performance of DMOC and its 

relationship to other existing optimal control methods are investigated. 
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1 Introduction 

In order to solve optimal control problems for mechanical systems, this paper links two important areas of re- 
search: optimal control and variational mechanics. The motivation for combining these fields of investigation 
is twofold. Besides the aim of preserving certain properties of the mechanical system for the approximated 
optimal solution, optimal control theory and variational mechanics have their common origin in the calculus 
of variations. In mechanics, the calculus of variations is also fundamental through the principle of stationary 
action also called Hamilton's principle. When applied to the action of a mechanical system, this principle 
yields the equations of motion for that system — the Euler- Lagrange equations. In optimal control theory 
the calculus of variations also plays a fundamental role. For example, it is used to derive optimality con- 
ditions via the Pontryagin maximum principle. In addition to its importance in continuous mechanics and 
control theory, the discrete calculus of variations and the corresponding discrete variational principles play 
an important role in constructing efficient numerical schemes for the simulation of mechanical systems and 
for optimizing dynamical systems. 

Discrete Optimal Control and Discrete Variational Mechanics. The theory of discrete variational 
mechanics has its roots in the optimal control literature of the 1960s; sec for example Cadzow [1970]; 
Hwang and Fan [1967]; Jordan and Polak [1964]. Specifically, Cadzow [1970] developed a discrete calculus 
of variations theory in the following way: A function is introduced which depends on a sequence of numbers, 
e.g. a sequence of times. A minimizing sequence necessarily satisfies a second-order difference equation, 
which is called the discrete Euler equation in reminiscence of its similarity with the Euler equation of the 
classical calculus of variations. 

An application of the discrete calculus of variations to an optimal control problem leads to a so called 
direct solution method. In this, one transforms the optimal control problem into a finite dimensional equality 
constrained nonlinear optimization problem via a finite dimensional parametrization of states and controls. 
In contrast, indirect methods (see §3.3 for an overview) are based on the explicit derivation of the necessary 
optimality conditions via the Pontryagin maximum principle. 

On the other hand, the theory of discrete variational mechanics describes a variational approach to 
discrete mechanics and mechanical integrators. The application of a discrete version of Hamilton's principle 
results in the discrete Euler-Lagrange equations. Analogous to the continuous case, near conservation of 
discrete energy, discrete momentum maps related to the discrete system's symmetries and the discrete sym- 
plectic form can be shown. This is due to the discretization of the variational structure of the mechanical 
system directly. Early work on discrete mechanics was often independently done by Cadzow [1973]; Lee 
[1983, 1987]; Logan [1973]; Macda [1980, 1981a,b]. In this work, the role of the discrete action sum, the 
discrete Euler-Lagrange equations and the discrete Noether's theorem were clearly understood. The varia- 
tional view of discrete mechanics and its numerical implementation is further developed in Wcndlandt and 
Marsden [1997a, b] and then extended in Bobcnko and Suris [1999a,b]; Kane, Marsden, and Ortiz [1999]; 
Kane, Marsden, Ortiz, and West [2000]; Marsden, Pckarsky, and Shkollcr [1999a,b]. The route of a vari- 
ational approach to symplectic-momentum integrators has been taken by Suris [1990], MacKay [1992]; see 
the review by Marsden and West [2001] and references therein. In this review a detailed derivation and 
investigation of these variational integrators for conservative as well as for forced and constrained systems 
is given. 
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Combining Optimal Control and Variational Mechanics. The present paper concerns the optimal 
control of dynamical systems whose behavior can be described by the Lagrange- d'Alembert principle. To 
numerically solve this kind of problem, we make use of the discrete calculus of variations only, that means 
we apply the discrete variational principle on two layers. On the one hand we use it for the description of 
the mechanical system under consideration, and on the other hand for the derivation of necessary optimality 
conditions for the optimal control problem. The application of discrete variational principles already on 
the dynamical level (namely the discretization of the Lagrange-d'Alembert principle) leads to structure- 
preserving time-stepping equations which serve as equality constraints for the resulting finite dimensional 
nonlinear optimization problem. The benefits of variational integrators are handed down to the optimal 
control context. For example, in the presence of symmetry groups in the continuous dynamical system, also 
along the discrete trajectory the change in momentum maps is consistent with the control forces. Choosing 
the objective function to represent the control effort, which has to be minimized is only meaningful if the 
system responds exactly according to the control forces. 

Related Work. A survey of different methods for the optimal control of dynamical systems described by 
ordinary differential equations is given in §3.3. However, to our knowledge, DMOC is the first approach 
to solutions of optimal control problems involving the concept of discrete mechanics to derive structure- 
preserving schemes for the resulting optimization algorithm. Since our first formulations and applications to 
space mission design and formation flying (Junge, Marsden, and Obcr-Blobaum [2005, 2006]; Junge and Ober- 
Blobaum [2005]; Obcr-Blobaum [2008]), DMOC has been applied for example to problems from robotics and 
biomechanics (Kanso and Marsden [2005]; Kobilarov, Desbrun, Marsden, and Sukhatmc [2007]; Kobilarov 
and Sukhatmc [2007]; Martin [2006]; Pckarck, Ames, and Marsden [2007]; Ross [2006]) and to image analysis 
(McLachlan and Marsland [2006]). From the theoretical point of view, considering the development of 
variational integrators, extensions of DMOC to mechanical systems with nonholonomic constraints or to 
systems with symmetries are quite natural and have already been analyzed in Kobilarov, Desbrun, Marsden, 
and Sukhatmc [2007]; Kobilarov and Sukhatmc [2007]. Further extensions are currently under investigation, 
for example DMOC for hybrid systems Pckarck, Ames, and Marsden [2007] and for constrained multi-body 
dynamics (see Leyendecker, Obcr-Blobaum, and Marsden [2008]; Leyendecker, Ober-Blobaum, Marsden, and 
Ortiz [2007]; Ober-Blobaum [2008]). DMOC related approaches are presented in Lee, McClamroch, and Leok 
[2006a,b]. The authors discretize the dynamics by a Lie group variational integrator. Rather than solving the 
resulting optimization problem numerically, they construct the discrete necessary optimality conditions via 
the discrete variational principle and solve the resulting discrete boundary value problem (the discrete state 
and adjoint system). The method is applied to the optimal control of a rigid body and to the computation 
of attitude maneuvers of a rigid spacecraft. 

Outline. In §2.1 and §2.2 we introduce the relevant concepts from classical variational mechanics and 
discrete variational mechanics following the work of Marsden and West [2001]. Especially, we focus on 
the Lagrangian and Hamiltonian description of control forces for the established framework of variational 
mechanics. Definitions and concepts of the variational principle, the Legendre transform, and Noether's 
theorem are readopted for the forced case in both the continuous as well as the discrete setting. In §3.1 and 
§3.2 we combine concepts from optimal control and discrete variational mechanics to build up a setting for 
the optimal control of a continuous and a discrete mechanical system, respectively. §2.3 describes the corre- 
spondence between the continuous and the discrete Lagrangian system as basis for a comparison between the 
continuous and the discrete optimal control problems in §3.3: We link both frameworks viewing the discrete 
problem as an approximation of the continuous one. The application of discrete variational principles for 
a discrete description of the dynamical system leads to structure-preserving time-stepping equations. Here, 
the special benefits of variational integrators are handed down to the optimal control context. These time- 
stepping equations serve as equality constraints for the resulting finite dimensional nonlinear optimization 
problem, therefore the described procedure can be categorized as a direct solution method. Furthermore, 
we show the equivalence of the discrete Lagrangian optimal control problems to those resulting from Runge- 
Kutta discretizations of the corresponding Hamiltonian system. This equivalence allows us to construct and 
compare the adjoint systems of the continuous and the discrete Lagrangian optimal control problem. In this 
way, one of our main results is related to the order of approximation of the adjoint system of the discrete 
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optimal control problem to that of the continuous one. With the help of this approximation result, we show 
that the solution of the discrete Lagrangian optimal control system converges to the continuous solution of 
the original optimal control problem. The proof strategy is based on existing convergence results of optimal 
control problems discretized via Runge-Kutta methods (Dontchev, Hager, and Veliov [2000]; Hagcr [2000]). 
§4.1 gives a detailed description of implementation issues of our method. Furthermore, in §4.2 we numerically 
verify the preservation and convergence properties of DMOC and the benefits of using DMOC compared to 
other standard methods to the solution of optimal control problems. 



2 Mechanical systems with forcing and control 
2.1 Variational Mechanics 

Our aim is to optimally control Lagrangian and Hamiltonian systems. For the description of their dynamics, 
we introduce a variational framework including external forcing resulting from dissipation, friction, loading 
and in particular control forces. To this end, we extend the notions in Marsden and West [2001] to Lagrangian 
control forces. 



Forced Lagrangian Systems. Consider an n-dimensional configuration manifold Q with local coordinates 
q = (q , . . . , q n ), the associated state space given by the tangent bundle TQ and a Lagrangian L : TQ — ► ML 
Given a time interval [0,T], we define the path space 

C(Q) = C([0, T], Q) = {q : [0, T] Q \ q is a C 2 curve} 

and the action map (3 : C — > K 

<S(q)= f L(q(t),q(t))dt. (2.1) 
Jo 

The tangent space T q C(Q) to C(Q) at the point q is the set of C 2 maps v q : [0, T] — + TQ such that tqov q = q, 
where tq : TQ — > Q is the natural projection. 

To define control forces for Lagrangian systems, we introduce a control manifold U C W 11 and define 
the control path space 

C(U) =C([Q,T],U) = {u : [0,T] -> U\u £ L°°}, 

with u(t) £ U also called the control parameter. Here, L°° denotes the space of essentially bounded, 
measurable functions equipped with the essential supremum norm. With this notation we define a Lagrangian 
control force as a map /l : TQ x U — > T*Q, which is given in coordinates as 

h ■ (q,<i,u) >-> (q,fL(q,q,u)), 

where we assume that the control forces can also include configuration and velocity dependent forces resulting 
e.g. from dissipation and friction. We interpret a Lagrangian control force as a family of Lagrangian forces 
with fixed curves u, that are fiber-preserving maps f^ '■ TQ —tT*Q over the identity id,Q, which we write in 
coordinates as 

Whenever we denote fL(q,q,u) as a one-form on TQ, we mean the family of horizontal one-forms f^(q,q) 
on TQ induced by the family of fiber-preserving maps /£. 

Given a control path u S C(U), the Lagrange- d'Alembert principle seeks curves q € C(Q) satisfying 

L(q(t),q(t))dt+ [ f L (q(t),q(t)Mt))-Sq(t)dt = Q, (2.2) 



where S represents variations vanishing at the endpoints. The second integral in (2.2) is the virtual work 
acting on the mechanical system via the force /l- Integration by parts shows that this is equivalent to the 
forced Euler- Lagrange equations 

^(Q,q) - ^ (|f + h(q,q,u) = 0. (2.3) 
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These equations implicitly define a family of forced Lagrangian flows and forced Lagrangian vector fields 
where for a fixed curve u € C(U), F£ : TQ x R — > TQ denotes the forced Lagrangian flow of the forced 
Lagrangian vector field X£ : TQ T(TQ). 

The one-form 0^ on TQ given in coordinates by 

FIT 

= of* <M> 

is called the Lagrangian one-form and the Lagrangian symplectic form SIl = d0£, is given in coordinates 
by 

Recall that the absence of forces, the Lagrangian symplectic form is preserved under the Lagrangian 
flow as (F^)* (£l L ) = fl Lj where F f L : TQ — > TQ is the map Fl at some fixed time t. 

Forced Hamiltonian Systems. Consider an n-dimensional configuration manifold Q, and define the 
phase space to be the cotangent bundle T*Q. The Hamiltonian is a function H : T*Q — > M. We will 
take local coordinates on T*Q to be {q,p) with q — (q , . . . , g") and p = (f>i, . . . ,p n ). Define the canonical 
one-form 6 on T*Q by 

@(Pq) ■ U Pq = (Pv T KQ ■ U P q )' ( 2 - 5 ) 

where p q G T*Q, u Pq G T Pq (T*Q), tt q : T*Q Q is the canonical projection, : T(T*Q) -> TQ is the 
tangent map of 7Tq and (•, •) denotes the natural pairing between vectors and covectors. In coordinates, we 
have 0(q,p) = Pidq\ The canonical two-form fl on T*Q is defined to be 

Q = -d6, 

which has the coordinate expression fl(q,p) — dq l A dpi. 

A Hamiltonian control force is a map : T* Q x J7 — > T* Q identified by a family of fiber preserving 
maps fjj :T*Q T*Q over the identity. Given such a control force, we define the corresponding family of 
horizontal one- forms f' H on T*Q by 

(/if)'(Pg) ' w p, = (fH(Pq): T ^ Q ■ w Pq ) 

for fixed curves u G C(f) where ttq : T*Q — * Q is the projection. This expression is reminiscent of definition 
(2.5) of the canonical one-form 9 on T*Q, and in coordinates it reads (/jj)'(g,p) ■ (Sq,Sp) = fjj{q,p) ■ Sq, 
thus the one-form is clearly horizontal. The forced Hamiltonian vector field Xjj for a fixed curve u is now 
defined by the equation 

and in coordinates this gives the well-known forced Hamilton's equations 



x q(i,p) = -jjrfo.p), ( 2 - 6a ) 



^(ff,p) = -g-(«,p) + /H(«.p), (2.6b) 

which are the standard Hamilton's equations in coordinates with the forcing term added to the momentum 
equation. This defines the forced Hamiltonian flow Fjj : T*Q x M — > T*Q of the forced Hamiltonian vector 
field = (X^ X^) for a fixed curve u G C(E7). 

The Legendre Transform with Forces. Given a Lagrangian L, we can take the standard Legendre 
transform FL : TQ —fT*Q defined by 



¥L(v q ) ■ w q = — 



L{v q + ew g ), 

e=0 



2.1 Variational Mechanics 



6 



where v q ,w q £ T q Q, and which has coordinate form 

FL : (q,q) h-> (q,p) = iq, -qT^A) 
and relate Hamiltonian and Langrangian control forces by 

ft = r H o fl. 

If we also have a Hamiltonian H related to L by the Legendre transform as 

H(q,p)=FL(q,q)-q-L(q,q), 

then the forced Euler-Lagrange equations and the forced Hamilton's equations are equivalent. That is, if 
and Xfl are the forced Lagrangian and Hamiltonian vector fields, respectively, then (FL)*(X^) = X^. To 
see this, we compute 



dH dq dL dL dq 

-dq {q > P) = P -dq-dq iq > q) ~dq {q > q) dq 
dL 

= -dq^ q) 

= -p + f%oFL{q,q), 

= -p + f%(q,p), 
dH . dq dL dq 

-d P ^ p) = q+p -d p -M [q ' q) d P 



= q, 

where p = FL(q, q) defines q as a function of (q,p). 



Noether's Theorem with Forcing. A key property of Lagrangian flows is their behavior with respect to 
group actions. Assume a Lie group G with Lie algebra g acts on Q by the (left or right) action <fi : GxQ — > Q. 
Consider the tangent lift of this action to (j> T< ^ : G x TQ given by 4>^(v q ) = T(4> g ) ■ v q , which is 

4> TQ (9,(q,q))= ^(g,q),^(g,q)q^ . 
For £ S g the infinitesimal generators £q : Q — > TQ and £tq ■ TQ — > T(TQ) are defined by 

= ^ (Ml)) ■ Z, = ~ (cbJ Q (v q )) • 

and the Lagrangian momentum map Jl : TQ — > q* is defined to be 

Jl(« 9 ) ■ £ = ©L • £,TQ(v q ). 

If the Lagrangian is invariant under the lift of the action, that is we have L o cfy^Q = L for all g G G (we 
also say, the group action is a symmetry of the Lagrangian), the Lagrangian momentum map is preserved of 
the Lagrangian flow in the absence of external forces. We now consider the effect of forcing on the evolution 
of momentum maps that arise from symmetries of the Lagrangian. In Marsden and West [2001] it is shown 
that the evolution of the momentum map from time to time T is given by the relation 

■£= f f U L {q{t)A{t))-t Q {q{t))dt. (2.7) 
JO 

Equation (2.7) shows, that forcing will generally alter the momentum map. However, in the special case 
that the forcing is orthogonal to the group action, the above relation shows that Noether's theorem will still 
hold. 



(j L o(Fi:) T ) (g(0),q(0))-J L (9(0),g(0)) 
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Theorem 2.1 (Forced Noether's theorem). Consider a Lagrangian system L : TQ — > R with control forcing 
fh '■ TQ xU — > T*Q and a symmetry action <f> : GxQ — > Q swc/i i/ioi (f^(q, g), £g(g)) = /or aZZ (g, g) € TQ, 
it € C(Z7) and aZZ £ £ g. Then the Lagrangian momentum map Jl ■ TQ — > g* wwZZ &e preserved by the flow, 
such that 3 L o (T^)* = J L /or a// 1. 

2.2 Discrete Mechanics 

The Discrete Lagrangian. Again we consider a configuration manifold Q, and define the ("discrete") 
state space to be Q x Q. Rather than considering a configuration q and velocity q (or momentum p), we 
now consider two configurations go and gi , which should be thought of as two points on a curve q which are 
a time step h > apart, i.e. go ~ g(0) and gi w g(/i). The manifold Q x Q is locally isomorphic to TQ and 
thus contains the same amount of information. 

A discrete Lagrangian is a function L d : Q x Q — > R, which we think of as approximating the action 
integral along the exact solution curve segment g between g and gi : 

L d {q ,qi)~ f L{q{t),q{t))dt. (2.8) 
Jo 

We consider the grid {t k = kh \ k = 0, . . . , N}, Nh = T, and define the discrete path space 

C d {Q) = C d ({t k }» =0 ,Q) = {q d : {t k }to - Q}- 

We will identify a discrete trajectory q d e C d {Q) with its image g<j = {gfej^O' where q k = qd(tk)- The 
discrete action map <S d : C d (Q) — > M along this sequence is calculated by summing the discrete Lagrangian 
on each adjacent pair and defined by 

N-l 

©d(gd) = X] L d{qk,qk+i)- 

k=0 

As the discrete path space C d is isomorphic to Q x • • • x Q (N + 1 copies), it can be given a smooth product 
manifold structure. The discrete action (3 d inherits the smoothness of the discrete Lagrangian L d . 

The tangent space T qd C d (Q) to C d (Q) at q d is the set of maps v qd : {tk}k=o ~ * sucn that r q ov qd = q d , 
which we will denote by v qd = {(qk, v k)}k=cr 

To complete the discrete setting for forced mechanical systems, we present a discrete formulation of the 
control forces introduced in the previous section. Since the control path u : [0,T] — > U has no geometric 
interpretation, we have to find an appropriate discrete formulation to identify a discrete structure for the 
Lagrangian control force. 

Discrete Lagrangian Control Forces. Analogous to the replacement of the path space by a discrete path 
space, we replace the control path space by a discrete one. To this end we consider a refined grid At, generated 
via a set of control points < c\ < ■ ■ ■ < c s < 1 as At — {the = tk + ceh \ k = 0, . . . , N — 1, 1 = 1, . . . , s}. 
With this notation the discrete control path space is defined to be 

C d (U) = C d (At, U) = {u d : At -> U}. 

We define the intermediate control samples u k on *fc+i] as u k = (u k i, ■ ■ ■ ,u ks ) € U s to be the values 
of the control parameters guiding the system from q k = q d (tk) to g^+i = q d (tk+i), where Uki = u d (tki) for 

le{i,...,s}. 

With this definition of the discrete control path space, we take two discrete Lagrangian control forces 
fd'fd '■ Q x Q x U s — » T* Q , given in coordinates as 



fd(Qk,Qk+i,Uk) = {qk+i, fd{qk,Qk+i,Uk)), 
fd(Qk,Qk+i,Uk) = (Qk,fd(Qk,Qk+i,Uk)), 



(2.9) 
(2.10) 
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also called left and right discrete forces. 1 Analogously to the continuous case, we interpret the two discrete 
Lagrangian control forces as two families of discrete fiber-preserving Lagrangian forces /^ fc,± :QxQ^T*Q 
in the sense that ttq a f^ k ' = ttq with fixed u k € U s and with the projection operators 7Tq : Q x Q — > 

Q, (qk,qk+i) i— * <Zfc+i and -Kq : Q x Q — > Q, [q kl q k +i) Qk- We combine the two discrete control forces to 
give a single one- form f^ k : Q x Q — > T* (Q x Q) defined by 



fd k {Qk,qk+i) ■ (Sq k ,Sq k+1 ) = f^ k ' + (qk,qk+i) • %+i + (?fc,?fc+i) ■ <fyk, 



(2.11) 



where f d (q k ,qk+i,u k ) denotes the family of all one-forms f d k {q k ,q k +i) with fixed u k € U s . To simplify the 
notation we denote the left and right discrete forces by f k := f^(qk, Sfc+i, u k), respectively, and the pair 
consisting of both by f k := f d (qk, Qk+i, Uk)- 

Referring to Figure 2.1, we interpret 
the left discrete force f k _ 1 as the force re- 
sulting from the continuous control force 
acting during the time span [t k _i,t k ] on 
the configuration node q k . The right dis- 
crete force is the force acting on q k re- 
sulting from the continuous control force 
during the time span [t k ,t k +i]. 




Figure 2.1: Left and right discrete forces. 

The Discrete Lagrange-d'Alembert Principle. As with discrete Lagrangians, the discrete control 
forces also depend on the time step h, which is important when relating discrete and continuous mechanics. 
Given such forces, we modify the discrete Hamilton's principle, following Kane, Marsdcn, Ortiz, and West 
[2000], to the discrete Lagrange-d'Alembert principle, which seeks discrete curves {q k } k= o that satisfy 



N-l 



N-l 



S ^2 L d (qk,qk+i) + ^2 [f d (qk,qk+i,u k ) ■ Sq k + f^(q k ,q k+1 ,u k ) ■ 5q k+1 ] = 



(2.12) 



fe=0 



k=0 



for all variations {&qk} k= o vanishing at the endpoints. This is equivalent to the forced discrete Euler- Lagrange 
equations 

D 2 L d (qk-i,qk) + D 1 L d (q k ,q k+1 ) + f^{q k -i,qk,u k -i) + fj(qk,qk+i,u k ) = 0. (2.13) 

These equations implicitly define the forced discrete Lagrangian map p™*- 1 ' Uk ■ q x q _> q x q f or fixed 
controls u fe _i,« fe S U s , mapping (q k -i,q k ) to (q k ,q k+1 ). 

The discrete Lagrangian one-forms ej d and e>i d are in coordinates 



L d (9o,'7i) = D 2 L d (q ,q 1 )dq 1 = 
L d (<7o,<7i) = -D 1 L d (q ,q 1 )dq 



dL d 
dq\ 



dq\, 



dL d 
dql 



dq l . 



(2.14a) 
(2.14b) 



In the absence of external forces, the discrete Lagrangian maps inherit the properties of symplectic 
preservation from the continuous Lagrangian flows. That means the discrete Lagrangian symplectic form 



flL d — dOj d = d6 Ld , with coordinate expression 



^Ld(<7o,gi) = 



d 2 L d 



dql A dgj 



dq^dqi 

is preserved under the discrete Lagrangian map as (Fi, d )*(QL d ) = r^ d , if no external forcing is present. 



Observe that the discrete control force is now dependent on the discrete control path. 
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The Discrete Legendre Transforms with Forces. Although in the continuous case we used the stan- 
dard Legendre transform for systems with forcing, in the discrete case it is necessary to take the forced 
discrete Legendre transforms 

¥ f+ L d : (q ,q 1 ,u Q ) h-> (qi,pi) = {qi, D 2 L d (q ,q 1 ) + f d (qo, qi, Uq)), (2.15a) 
¥ f ~L d : (q 0) qi,u ) i-» (qo,Po) = (go, -DiL d (q 0) g x ) - f d {q , qi, u )). (2.15b) 

Again, we denote with ¥' ± L d t ° the forced discreteLegendre transforms for fixed controls uo <E U s . Using 
these definitions and the forced discrete Euler-Lagrange equations (2.13), we can see that the corresponding 
forced discrete Hamiltonian map 

is given by the map : (qo,Po) ^ (li,Pi), where 

Po = -DiL d {q ,qi) - f d °^(qn,qi), (2.16a) 
Pl = D 2 L d (q Q ,q 1 ) + f^ + (q ,q 1 ), (2.16b) 

which is the same as the standard discrete Hamiltonian map with the discrete forces added. 

Figure 2.2 shows that the following two definitions of the forced discrete Hamiltonian map 

puo = ¥ f± L ui puo.ux Qp/i^uo)-^ (2.17 a ) 

puo = F /+i«o (F/-L7)" 1 , (2.17b) 

are equivalent with coordinate expression (2.16). Thus from expression (2.17b) and Figure 2.2, it becomes 
clear, that the forced discrete Hamiltonian map that maps (qo,Po) to (qi,pi) depends on uq only. 




Figure 2.2: Correspondence between the forced discrete Lagrangian and the forced discrete Hamiltonian 
map. 



The Discrete Noether Theorem with Forcing. As in the unforced case, we can formulate a discrete 
version of the forced Noether 's theorem (for the derivation see for example Marsden and West [2001]). To 
this end, the discrete momentum map in presence of forcing is defined as 

J L + d (q ,qi)-Z=(V f+ L u d °( qQ , qi )^ Q ( qi )), 
j{-(qo,qi)-Z = (Ff-L u d °(q , qi ),Z Q (q )). 
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The evolution of the discrete momentum map is described by 

JV-l 

N-l 



JiMO - J L d \ (90,gi)-e=E /?(ft,ft + l)^QxQfe,9Hl)- (2-18) 

fc=0 

Again, in the case that the forcing is orthogonal to the group action we have the unique momentum map 
j{ d : Q x Q — > q* and it holds: 

Theorem 2.2 (Forced discrete Noether's theorem). Consider a discrete Lagrangian system Ld : Q x Q — > R 
with discrete control forces fd>fd '■ Q x Q x U s ^ T*Q and a symmetry action <f> : G x Q — > Q such that 
(fd k > Cqxq) = /or all£ e Q and £ U s , k £ {0, . . . , TV— 1}. TTien ifte discrete Lagrangian momentum map 
^L d •QxQ ~* &* w ^ be preserved by the discrete Lagrangian evolution map, such that J £ o p^ Uk + 1 = j£ 

2.3 The Discrete vs. the Continuous Lagrangian Systems 

In this section, we relate the continuous and the discrete Lagrangian system. First, along the lines of 
Marsden and West [2001], we define expressions for the discrete mechanical objects that exactly reflect the 
continuous ones. Based on the exact discrete expressions, we determine the order of consistency concerning 
the difference between the continuous and the discrete mechanical system. 

Exact Discrete Lagrangian and Forcing. Given a regular Lagrangian L : TQ — > K and a Lagrangian 
control force /l : TQ x U — ► T*Q, we define the exact discrete Lagrangian L^ : Q x Q x K — > K and the 
exact discrete control forces f^ + , /<f ~ : Q x Q x C([0, ft] , £/) xl-> T*Q to be 

L$(q ,q l} h)= / L(g(t),g(t))di, (2.19) 
Jo 

f? + (q ,qun Q ,h)= [ /L(9(t),g(t),«(t))-^di, (2.20) 

Jo "91 

/ d E -(g ,9i,u ,ft) = /" /z(g(*), ?(*),«(*)) -^dt, (2.21) 

Jo "90 

with Ufc € C([fcft, (fe+l)ft], C/) and where g : [0, ft] — > Q is the solution of the forced Euler-Lagrange equations 
(2.3) with control function u : [0, ft] — > U for L and /l satisfying the boundary conditions g(0) = go and 
q(h) = qi. Observe, that the exact discrete control forces depend on an entire control path in contrast to 
the continuous control forces. Consequently, the exact forced discrete Legendre transforms are given by 

F'+Lf (ft,, 9i, uo, ft) - (qi,D 2 L%(q , q u h) + /J+(g , q x , u 0) ft)), 
^ f ~ L d (9o, 9i, u , ft) = (g ,--Di£f (go,9i,ft) - /<f~(9o, 9i,u , ft))- 

As in §2.2 F/±Lf Ufc (9fe,9 fc+ i) and /f , gfe+i) denote the exact discrete forces and the exact forced 

discrete Legendre transforms for a fixed U& S C([fcft, (fe + l)ft], f ). The next lemma is based on a result in 
Marsden and West [2001] (Lemma 1.6.2) extended to the presence of control forces and establishes a special 
relationship between the Legendre transforms of a regular Lagrangian and its corresponding exact discrete 
Lagrangian. This also proves that exact discrete Lagrangians are automatically regular. 

Lemma 2.3. A regular Lagrangian L and the corresponding exact discrete Lagrangian L® have Legendre 
transforms related by 



F /+ Lf (g , 9i, u , ft) = FZ(g ,i(ft), g ,i(ft)), 
F'-Lf (g , 9i, u , ft) = FL(g ,i(0), go,i(0)), 



for sufficiently small ft and close go,9i G Q- Here go.i denotes the solution of the corresponding Euler- 
Lagrange equations with g(0) = go,g(ft) = gi- 
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Proof. Analogous to the proof in Marsden and West [2001] (Lemma 1.6.2) for the unforced case we begin 
with Ff~Lj and compute 



F/ L d(<10,qi,U-o,h) 



9L dg QA dL dq Q ,i 
dq dq dq dq 



dt 



h 



dq 



dt 



h r 



dL _ d_dL 
dq dt dq 



h 



dqo, 



0q { 



idt- 



dL dq ,i 



dq dq a 



using integration by parts. Since qo,i(t) is a solution of the forced Euler-Lagrange equations the first term 
is zero. With qo,i(0) = 1o anc ^ Qo,i(h) = 9i for the second term we get 

— — (0)=id and — — (h)=0. 
dq dq 

Substituting these into the above expression for ¥^~L^ now gives 

dL 

Wf-L$(q , qi ,u ,h) = -gr (5o,i (0), Qo,i (0)), 

which is simply the definition of FL(go,i(0), go,i(0)). 

The result for F^+Lf can be established by a similar computation. ■ 

Combining this result with the diagram in Figure 2.2 gives the commutative diagram shown in Figure 
2.3 for the exact discrete Lagrangian and forces. The diagram also clarifies the following observation, that 
was already proved in Marsden and West [2001] (Theorem 1.6.3) for unforced systems and can now be 
established for the forced case as well: Consider the pushforward of both, the continuous Lagrangian and 
forces and their exact discrete Lagrangian and discrete forces to T*Q, yielding a forced Hamiltonian system 
with Hamiltonian H and a forced discrete Hamiltonian map , respectively. Then, for a sufficiently 
small time step the forced Hamiltonian flow map equals the pushforward discrete Lagrangian map: 



{F u 0) h 



F u 



Order of Consistency. In the previous paragraph we observed that the exact discrete Lagrangian and 
forces generate a forced discrete Hamiltonian map that exactly equals the forced Hamiltonian flow of the 
continuous system. Since we are interested in using discrete mechanics to reformulate optimal control 
problems, we generally do not assume that Ld and L or H are related by (2.19). Moreover, the exact 
discrete Lagrangian and exact discrete forces are generally not computable. In this section we determine the 
error we obtain by using discrete approximations for the Lagrangian and the control forces. 

Rather than considering how closely the trajectory of F matches the exact trajectory given by Fjj, as 
is done in forward error analysis, we use the concept of variational error analysis as introduced in Marsden 
and West [2001]. In this context we consider how closely a discrete Lagrangian matches the exact discrete 
Lagrangian given by the action. For forced systems, we additionally take into account how closely the 
discrete forces match the exact discrete forces. As was stated in the last section, if the discrete Lagrangian 
is equal to the action and the discrete forces are given by (2.20) and (2.21), then the corresponding forced 
discrete Hamiltonian map F^ will exactly equal the forced flow F^ k . Usually, this is just an approximation, 
therefore we define the local variational error as follows (see Marsden and West [2001]). For the following 
we assume Q to be a normed vector space equipped with the norm || • ||. 

Definition 2.4 (Order of consistency). A given discrete Lagrangian Ld is of order r, if for a fixed curve 
u G C(U) there exist an open subset V v C TQ with compact closure and constants C v > and h v > such 
that 

\L d (q(0),q(h), h) - L^{q{Q),q{h),h)\ < C v h r+1 (2.22) 

for all solutions q{t) of the forced Euler-Lagrange equations with initial condition (q°,q°) G V v and for all 
h < h v . Analogously, a given discrete force f^ ' is of order r, if there exist an open subset V w C TQ with 
compact closure and constants C w > and h w > such that 

||/ d uo,± (3(0), gW,^)-/^ ^^),^),/!)!! < C w h r+1 (2.23) 
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(<7o,9i) h 




(90, Po) h 



(9o,9o) h 



^uo = (j pu 0) ^ 



{F u 0) h 



-* (9i, Pi) h 



-> (9i, 9i) h 



{F u 1)h 



~> (92,P2) 



F, 



-> (92,92) 



Figure 2.3: Correspondence between the exact discrete Lagrangian and forces and the continuous forced 
Hamiltonian flow. 



for all solutions q(t) of the forced Euler-Lagrange equations with initial condition (9,9) € V w for all 
h < h w and with fixed uq G U s and uq G C([0,h],U). The discrete Legendre transforms ¥ + L d ° and ¥~ L u d 
of a discrete Lagrangian are of order r if there exist an open subset Vt C TQ with compact closure and 
constants Ct > and hf > such that 

\\¥ ± L^( q (0), q (h) 1 h)-¥ ± L^( q (0), q (h),h)\\ < C f h r+1 (2.24) 

for all solutions q(t) of the forced Euler-Lagrange equations with initial condition (9,9) G Vf and for all 
h < hf and with fixed uq G U s and vlq G C([0, h],U). 

To give a relationship between the orders of a discrete Lagrangian, discrete forces, the forced discrete 
Legendre transforms, and their forced discrete Hamiltonian maps, we first have to introduce what we mean 
by equivalence of discrete Lagrangians: L l d is equivalent to L 2 d if their discrete Hamiltonian maps are equal. 
For the forced case, we say analogously, that for fixed Uk G II s the discrete pair (L d , f d k ' 1 ) is equivalent to 
the discrete pair (L d , f d k ' ) if their forced discrete Hamiltonian maps are equal, such that F"f = . With 

pu, = F /+ £ «»,i (fZ-LJ*' 1 ) -1 , it follows that if (L d , f d kA ) and {L% /J fc ' 2 ) are equivalent, then their forced 
discrete Legendre transforms are equal. Thus, equivalent pairs of discrete Lagrangians and control forces 
generate the same integrators. 

The following theorem is an extended version of the unforced case stated in Marsdcn and West [2001] 
(Theorem 2.3.1). 

Theorem 2.5. Given a regular Lagrangian L, a Lagrangian control force fi and a corresponding Hamilto- 
nian H with Hamiltonian control force fn- Then for a fixed curve u G C(U) and a corresponding fixed control 
sequence Uk G II s the following statements are equivalent for a discrete Lagrangian Ld and the discrete forces 

fh- 

(i) the forced discrete Hamiltonian map for (Ld, /j*' ) is of order r, 
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(ii) the forced discrete Legendre transforms of (Ld, /^ fc ' ± ) are of order r, 

(Hi) (Ld, fd h ' ) * s equivalent to a pair of discrete Lagrangian and discrete forces, both of order r. 

Proof. Considering the proof of the unforced version in Marsden and West [2001] (Theorem 2.3.1), this 
proof is straightforward by taking into account the discrete Lagrangian control forces. For an alternative 
approach of the proof of Theorem 2.3.1 in Marsden and West [2001] see Patrick and Cucll [2008]. ■ 

Note that, given a discrete Lagrangian and discrete forces, their order can be calculated by expanding 
the expressions for Ld(q(0),q(h),h) and fj 1 "^ in a Taylor series in h and comparing these to the same 
expansions for the exact Lagrangian and the exact forces, respectively. If the series agree up to r terms, then 
the discrete Lagrangian is of order r. Analogously, the discrete forces are of order r, if for both expansions 
the first r terms are identical. 



3 Optimal Control of a Mechanical System 
3.1 The continuous setting 

On the configuration space Q we consider a mechanical system described by a regular Lagrangian L : 
TQ — > M. Additionally, assume that a Lagrangian control force acts on the system and is defined by a map 
fh ■ TQ x U — > T*Q with : (q,q,u) i— > (q, /l(<?, q, it)) and u : [0,T] — > U, the time-dependent control 
parameter. Note that the Lagrangian control force may include both dissipative forces within the mechanical 
system and external control forces resulting from actuators steering the system. 



The Lagrangian Optimal Control Problem. We now consider the following optimal control problem: 
During the time interval [0,T], the mechanical system described by the Lagrangian L is to be moved on a 
curve q from an initial state (q(0),q(0)) = (q°,q°) € TQ to a final state. The motion is influenced via a 
Lagrangian control force /& with control parameter u such that a given objective functional 

J(q,u) = ( T C(q(t),q(t),u(t))dt + $(q(T),q(T)) (3.1) 



is minimized. Here C : TQ X U — ► K and <i> : TQ — > E (Mayer term) are continuously differentiablc 
cost functions. The final state (q(T),q(T)) is required to fulfil a constraint r(q(T),q(T),q T ,q T ) = with 
r : TQ x TQ ^ R n - and (g T , q T ) € TQ given. 

The motion of the system is to satisfy the Lagrange-d'Alembert principle, which requires that 

L(q(t),q(t))dt+ [ f L (q(t),q(t),u(t))-Sq(t)dt = (3.2) 



JO 



for all variations Sq with Sq(0) = Sq(T) = 0. In many cases, one encounters additional constraints on the 
states and controls given by h(q(t) , q(t) , u(t)) > with h : TQ x U 
To summarize, we are faced with the following 

Problem 3.1 (Lagrangian optimal control problem (LOCP)). 



subject to 



min J(<7, it) (3.3a) 

qec(Q),uec(U) 



L(q(t),q(t))dt+ [ f L (q(t),q(t),u(t))-5q(t) dt = 0, (3.3b) 
JO 

(g(0),g(0)) = (gV), (3.3c) 

h(q(t),q(t),u(t)) > 0, te[0,T], (3.3d) 

r(q(T),q(T),q T ,q T ) = 0. (3.3e) 



3.1 The continuous setting 



14 



The interval length T may cither be fixed, or appear as degree of freedom in the optimization problem. 

Definition 3.2. A curve (q,u) G C(Q) xC(U) is admissible (or feasible ), if it fulfills the constraints (3.3c)- 
(3.3e). The set of all admissible (feasible) curves is the admissible (feasible) set of Problem 3.1. An admissible 
(feasible) curve (q*,u*) is an optimal solution of Problem 3.1, if 

J(q*,u*)<J(q,u) (3.4) 

for all admissible (feasible) curves (q,u). An admissible (feasible) curve (q*,u*) is a local optimal solution. 
if (3.4) is fulfilled in a neighborhood of (q*,u*). The function q* is called (locally) optimal trajectory, and 
u* is the (locally) optimal control. 



The Hamiltonian Optimal Control Problem. We now formulate the problem using the Hamiltonian 
variant for the system dynamics. This is equivalent to the Lagrangian formulation as we have seen in 
paragraph §2.1 on the Legendre transform with forces. 

For a set 1Z — {(q, q) G TQ \ g(q, q) > 0} determined via a constraint g : TQ — > W lg on TQ we obtain 
the corresponding set in the cotangent bundle as 1Z — {(q,p) G T*Q | g(q,p) > 0} with g = g o (FL)^ 1 . 
Correspondingly, the optimal control problem in the Hamiltonian formulation reads as follows: 

Problem 3.3 (Hamiltonian optimal control problem (HOCP)). 

J(q,p,u)= [ C(q(t),p(t),u(t))dt + $(q(T),p(T)) (3.5a) 



(q,p)€T'C(Q),u£C(U) 



subject to 



q(t) = W p H(q(t),p(t)), (3.5b) 

p(t) = -V q H(q(t),p(t)) + f H (q(t),p(t),u(t)), (3.5c) 

(?(0),p(0)) = (q°,p°), (3.5d) 

h(q(t),p(t),u(t)) > 0, t G [Q, T] (3.5e) 

f(q(T),p(T),q T ,p T ) = 0. (3.5f) 

where (q T ,p T ) = ¥L(q T ,q T ), p(0) = D 2 L(q(0),q(0)), p° = D 2 L(q°,q ) and <l = $0 (FL)" 1 etc. 



Necessary Optimality Conditions. In this paragraph, we introduce necessary conditions for the opti- 
mality of a solution (q*,u*) G C(Q) x C(U) to Problems 3.1 and 3.3. We restrict ourselves to the case of 
problems with the controls pointwise constrained to the (nonempty) set (/={ii£ M. n " \ h{u) > 0} and fixed 
final time T. With x — (q,p) and 

f(x,u) = (W p H(q,p),-W q H(q,p) + f H (q,p,u)) 

we can rewrite (3.5) as 

min J(x,u)= [ C(x(t),u(t))dt+$(x(T)) (3.6a) 

subject to 

x{t) = f(x(t),u(t)), (3.6b) 

x(0) — x , (3.6c) 

u(t) G U, te[0,T] (3.6d) 

f{x{T),x T ) = 0. (3.6e) 

Necessary conditions for optimality of solution trajectories ?/(•) = {x{-),u{-)) can be derived based on varia- 
tions of an augmented cost function, the Lagrangian of the system: 

£(77, A) = [ T C(x(t) lU (t)) + \ T (t)-(x-f(x(t),u(t))) dt + Q(x(T)), (3.7) 
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where the variables A^, i — 1, . . . ,n x , are the adjoint variables or the costates. A point (n* , A*) is a saddle 
point of (3.7), if for all r] and A it holds 

£{V, A*) < C{rf, A*) < C(rf, A). 

The function 

u, A) := -C(a;, it) + A T • /(x, u) (3.8) 

is called the Hamiltonian of the optimal control problem. When setting variations of C with respect to r\ 
and A to zero, the resulting Euler-Lagrange equations provide necessary optimality condition for the optimal 
control problem (3.6). Formally, one obtains the following celebrated Theorem (cf. Pontryagin, Boltyanski, 
Gamkrelidze, and Misccnko [1962]): 

Theorem 3.4 (Pontryagin Maximum Principle). Let (x*,u*) be an optimal solution to (3.6). Then there 
exists a piecewise continuously differentiate junction A : [0, T] — > R" x and a vector a G M™ r such that 

H{x*(t),u*(t),X(t)) = maxH(x(t),u,\(t)) ie[0,T], (3.9a) 

and A solves the following initial value problem: 

A(T) = V x *(x*(T))- V x r(x*{T),x T )a, (3.9b) 
A = -V x H(x*,u*,\). (3.9c) 

Transformation to Mayer Form. Within the previous sections we introduced optimal control problems 
in Bolza form, in which the objective functional consists of a cost functional of integral form and a final 
point constraint. For the error analysis in §3.3 it will be useful to transform the problem into Mayer form, 
in which the objective functional consists of the final point constraint only. To this end we introduce a new 
state variable as 

z(t) := f C(q(T),q(r),u(T))dT, < t < T, (3.10) 
Jo 

resp. y(t) := C (q(r) , p(r) , u(t)) dr, < t < T for the Hamiltonian form. By extension of the state space 
from TQ to TQ x R, and from T*Q to T*Q x K, respectively, the new objective function in Mayer form 
reads 

J(q,q,u) = z(T) + $(q(T),q(T)) 

resp. J(q,p,u) = y(T) + <&(q(T),p(T)). Equation (3.10) is typically adjoint to the problem description as an 
additional differential equation of the form 

z = C{q,q,u), Z(0) = 0, (3.11) 

resp. ij = C(q,p, u), y(0) = 0. 
3.2 The discrete setting 

For the numerical solution we need a discretized version of Problem 3.1. To this end we formulate an optimal 
control problem for the discrete mechanical system described by discrete variational mechanics introduced 
in §2.2. In §3.3 we show how the optimal control problem for the continuous and the discrete mechanical 
system are related. 

To obtain a discrete formulation, we replace each expression in (3.3) by its discrete counterpart in terms 
of discrete variational mechanics. As described in §2.2, we replace the state space TQ of the system by Q x Q 
and a path q : [0, T] — > Q by a discrete path q& : {0, h, 2h, . . . , Nh = T} — > Q with q^ = q^(kh). Analogously, 
the continuous control path u : [0,T] — > U is replaced by a discrete control path : At — » U (writing 
u k = (u d (kh + c t h)) s e=1 e U s ). 



3.2 The discrete setting 



16 



The Discrete Lagrange-d'Alembert Principle. Based on this discretization, the action integral in 
(3.2) is approximated on a time slice [kh, (fc + l)h] by the discrete Lagrangian L d : Q x Q — ► R, as in 
equation (2.8), 

L d (q k ,q k+1 ) w / L(q(t),q(t))dt, 

Jkh 

and likewise the virtual work by the left and right discrete forces, 

Ak+l)h 

fk ■ + fk ■ tlk+l « / f L (q(t),q(t),u(t)) ■ 6q(t) dt, 

Jkh 

where /-,/+€ T*Q. 

As introduced in equation (2.12), the discrete version of the Lagrange-d'Alembert principle (3.2) requires 
one to find discrete paths {q k } k=0 sucn that for all variations {Sq k }^ =0 with Sqo = 5qw = 0, one has 

JV-l 2V-1 
S L d (q k ,q k+1 ) + J2 [fk ■ Sqk + ft • Sq k+1 ] = 0, (3.12) 
fc=0 fc=0 

or, equivalently, the forced discrete Euler-Lagrange equations 

DiL A (q k -i, q k ) + D x L d (q k , q k+1 ) + f+_ t + f~ = 0, (3.13) 

k= 1,...,N- 1. 

Boundary Conditions. In the next step, we need to incorporate the boundary conditions q(0) = q°,q(0) = 
q° and r(q(T),q(T),q T ,q T ) = into the discrete description. Those on the configuration level can be used 
as constraints in a straightforward way as qo = q . However, since in the present formulation velocities are 
approximated in a time interval [t k ,t k+ i] (as opposed to an approximation at the time nodes), the velocity 
conditions have to be transformed to conditions on the conjugate momenta. These are defined at each time 
node using the discrete Legendre transform. The presence of forces at the time nodes has to be incorporated 
into that transformation leading to the forced discrete Legendre transforms L d and V^ + L d defined in 
(2.15). Using the standard Legendre transform FL : TQ — > T*Q, (q, q) i— > (q,p) = (q, D2L(q, q)) leads to the 
discrete initial constraint on the conjugate momentum 

D 2 L(q°, q°) + D^fa, qi ) + fj(q , qi,u ) = 0. 

As shown in the previous section, we can transform the boundary condition from a formulation with con- 
figuration and velocity to a formulation with configuration and conjugate momentum. Thus, instead of 
considering a discrete version of the final time constraint r on TQ we use a discrete version of the final time 
constraint f on T*Q. We define the discrete boundary condition on the configuration level to be 

r d : Q x Q x U s x TQ -> R n - , 

r d {qN~i,qN,u N -i 1 q T ,q T ) = f {¥ !+ L d (q N ^ l7 q N , ujv-i), JFL(<? T , q T )) , 

with (q Nl p N ) = F /+ i d (g A r_i,gr A r,M W -i) and (q T ,p T ) = FL(q T ,q T ), that is p N = D 2 L d (q N ^ 1 ,q N ) + 
fctiQN-iiqNiUN-x) and p T = D 2 L(q T ,q T ). 

Notice that for the simple final velocity constraint 

r(q(T) iq (T),q T ,q T ) = q(T)-q T , 

we obtain for the transformed condition on the momentum level r(q(T),p(T), q T ,p T ) = p(T)—p T the discrete 
constraint 

- D 2 L(q T ,q T ) + D 2 L d (q N - 1 , q N ) + f£(q N -i,q N ,u N -i) = 0. (3-14) 
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Discrete Path Constraints. Opposed to the final time constraint we approximate the path constraint 
in (3.3d) on each time interval [i/c,ifc+i] rather than at each time node. Thus, we maintain the formulation 
on the velocity level and replace the continuous path constraint h{q{i) , q(t) , u(t)) > by a discrete path 
constraint h d : Q x Q x II s — > R srih with 

h d (q k ,q k+1 ,u k ) > 0, fc = 0, ...,N-l. 



Discrete Objective Function. Similar to the Lagrangian we approximate the objective functional in 
(3.1) on the time slice [kh, (k + l)h] by 

Ak + l)h 

C d {q k ,q k +i,u k ) ~ / C(q(t),q(t),u(t))dt. 

Jkh 

Analogously to the final time constraint, we approximate the final condition via a discrete version <5>d : 
Q x Q x U" — *■ K yielding the discrete objective function 

AT-l 

Jd(qd,Ud) = ^2 c d(qk,qk+i,u k ) +®d{<iN-i,qN,UN-i)- 

k=0 

The Discrete Optimal Control Problem. In summary, after performing the above discretization steps, 
one is faced with the following discrete optimal control problem. 

Problem 3.5 (Discrete Lagrangian optimal control problem). 

min J d {q d ,u d ) (3.15a) 
qdM d 

subject to 

qo = q°, (3.15b) 

D 2 L(q°, q°) + D 1 L d (q , q x ) + / " = 0, (3.15c) 

D 2 L d (qk-l,qk) + D 1 Ld(q k ,qk+l) + ft 1 + fk=0, k = 1, . . . , N - 1, (3.15d) 

hd(qk,qk+i,Uk) > 0, ife = 0,...,JV-l, (3.15e) 

r d (qN-i,qN,u N - 1 ,q T ,q T ) = 0. (3.15f) 



Recall that the are dependent on u k € U s . To incorporate a free final time T as in the continuous 
setting, the step size h appears as a degree of freedom within the optimization problem. However, in the 
following formulations and considerations we restrict ourselves to the case of fixed final time T and thus 
fixed step size h. 



Necessary Optimality Conditions. The system (3.15) results in a constrained nonlinear optimization 
problem, also called a nonlinear programming problem, that is an objective function has to be minimized 
subject to algebraic equality and inequality constraints. Let £ be the set of parameters introduced by the 
discretization of an infinite dimensional optimal control problem. Then, the nonlinear programming problem 
to be solved is 

mine <p(£) 

* (3.16) 
subject to a(£) = 0, b(£) > 0, 

where (p : R" -> R, a : R" -> R m , and b : R n -> MP arc continuously differentiable. 

We briefly summarize some terminology: A feasible point is a point £ £ R" that satisfies a(£) = and 
&(0 > 0. A local minimum of (3.16) is a feasible point £* which has the property that v?(£*) < f° r an 
feasible points £ in a neighborhood of £*. A strict local minimum satisfies </?(£*) < for all neighboring 
feasible points £ 7^ £*. ^4cii?;e inequality constraints 6 act (£) at a feasible point £ are those components bj(£) 
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of &(£) with = 0. We subsume the equality constraints and the active inequalities at a point £ (known 
as active set) in a combined vector function of active constraints as 



«(0 

6 act (£) 



Note that the active set may be different at different feasible points. Regular points are feasible points £ that 
satisfy the condition that the Jacobian of the active constraints, Va(£) T , has full rank, that means that all 
rows of Va(£) T are linearly independent. 

To investigate local optimality in the presence of constraints, we introduce the Lagrangian multiplier 
vectors A £ M m and fi € M. p , that are also called adjoint variables, and we define the Lagrangian function C 
by 

£(£, A, (j,) := - A T a(£) - M T 6(£). (3.17) 

We now state a variant of the Karush-Kuhn- Tucker necessary conditions for local optimality of a point £*. 
These conditions have been derived first by Karush in 1939 (Karush [1939]) and independently by Kuhn and 
Tucker in 1951 (Kuhn and Tucker [1951]). For brevity, we restrict our attention to regular points only. 

Theorem 3.6 (Karush-Kuhn- Tucker conditions (KKT)). If a regular point £* £ K™ is a local optimum of 
the NLP problem (3.16), then there exist unique Lagrange multiplier vectors A* € K m and /i* € M. p such that 
the triple (£*,A*,/i*) satisfies the following necessary conditions: 

Vs£(f ,A*,M*)=0, 

«(D = o, 
b(C) > o, 
fi* > o, 
Mj^i(D = 0, 3 = 1,..., p. 

Proof. Sec for example Jarre and Stoer [2004]. ■ 

A triple (£*,A*,^t*) that satisfies the Karush-Kuhn- Tucker conditions is called a Karush- Kuhn- Tucker 
point (KKT point). 

Transformation to Mayer Form. As in the continuous setting, the discrete optimal control problem 
(3.15) can be transformed into a problem in Mayer form, that is the objective function consists of a final 
condition only. The transformation is performed on the discrete level to keep the Lagrangian structure of 
the original problem. We introduce new variables 



z Q = 0, Zj = y]Cd(qk,gk+i,Uk), i=l,...,N, 



k=0 

and rephrase the discrete problem (3.15) into a problem of Mayer type: 

min z N + ®<i(q N -i,q N ,u N -i) (3.18a) 
qd-ud 

subject to 

q = q°, (3.18b) 

z = 0, (3.18c) 

D 2 L(q ,q ) + D 1 L d (q ,q 1 ) + fo=0, (3.18d) 

D 2 L d (q k - 1 ,q k )+D 1 L d (q k ,q k+1 )+f+_ 1 +f-=0, k = l,...,N-l, (3.18e) 

h d (q k ,q k+1 ,u k )>0, k = 0,...,N-l, (3.18f) 

r d (qN-i,<lN,UN-i,q T A T ) = 0, (3-18g) 

Zk+i - z k - C d {q k ,q k+ i,u k ) = 0, k = 0, . . . , N - 1. (3.18h) 
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Thus equations (3.18c) and (3.18h) provide the corresponding discretization for the additional equation of 
motion (3.11) resulting from the Mayer transformation of the Lagrangian optimal control problem on the 
continuous level. 



Fixed Boundary Conditions. Consider the special case of a problem with fixed initial and final config- 
uration and velocities and without path constraints: 

min / C(q(t),q(t),u(t))dt (3.19a) 
q eC(Q),uec(U) J Q 

subject to 

Sf L(q(t),q(t))dt+ [ f L {q(t),q{t),u{t))-8q{t)dt = Q, (3.19b) 
Jo Jo 

g(0) = g°, g(0) = <?°, g(T) = g T , g(T) = q T . (3.19c) 

A straightforward way to derive initial and final constraints for the conjugate momenta rather than for 
the velocities from the variational principle directly is stated in the following proposition: 

Proposition 3.7. With (q°,p°) = ¥L{q°,q a ) and (q T ,p T ) = ¥L(q N ,q N ) equations (3.19b) and (3.19c) are 
equivalent to the following principle with free initial and final variation and with augmented Lagrangian 

L(q(t),q(t)) dt + p°(q(0)-q Q )-p T (q(T)-q T ))+ [ f L (q(t), q(t) , u(t)) ■ Sq(t) dt = 0. (3.20) 



Proof. Variations of (3.19b) with respect to q and zero initial and final variation 8q(0) = Sq(T) = together 
with (3.19c) yield 

±-^L(q(t),q(t)) - ^L(q(t),q(t)) = f L (q(t),q(t),u(t)), (3.21a) 

g(0) = g°, g(0) = g°, g(T) = q T ', g(T) = q T . (3.21b) 

On the other hand variations of (3.20) with respect to q and A = (p°,p T ) with free initial and final variation 
lead to 

^L(q(t),q(t)) - §- q L(q(t),q(t)) = f L (q(t),q(t),u(t)), (3.22a) 



§rL(q(t),q(t)) 
§.L(q(t),q(t)) 



= p T , (3.22b) 



= P°, (3.22c) 



t=T 



g(0) = g°, q(T) = q T . (3.22d) 

The Legendre transform applied to the velocity boundary equations in (3.21b) gives the corresponding 
momenta boundary equations (3.22b) and (3.22c). ■ 

On the discrete level we derive the optimal control problem for fixed initial and final configurations and 
velocities in an equivalent way. Thus, we consider the discrete principle with discrete augmented Lagrangian 

(N-i \ N-l 

Ld(qk,q k+ l)+P°(qo - q°)-p T (qN ~ q T )) + ]T [/,7 • 5q k + /+ • 5q k+1 ] = 0, (3.23) 
fe=0 / k=0 

which, with free initial and final variation Sqg and Sq^, respectively, is equivalent to 

JV-l AT-l 

S L d (q k ,q k+ i) + ttk ■ SQk + ft ■ <W] - 0, (3.24a) 

k=0 k=0 
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go = q°, P Q + D x L d {q Q , qi ) + / " = 0, (3.24b) 

q N = q T , - P T + D 2 L d (q N _ u q N ) + = 0, (3.24c) 

where the second equations in (3.24b) and (3.24c) are exactly the discrete initial and final velocity constraints 
derived in the remark containing equation (3.14) with p° = D2L(q°, q°) and p T = Z?2^(<7 T , q T )- 

We note that this derivation of the discrete initial and final conditions directly gives the same formulation 
that we found before by first transforming the boundary condition on the momentum level T*Q and then 
formulating the corresponding discrete constraints on Q x Q x U s . 

3.3 The Discrete vs the Continuous Problem 

This section gives an interpretation of the discrete problem as an approximation to the continuous one. In 
addition, we identify certain structural properties that the discrete problem inherits from the continuous 
one. We determine the consistency order of the discrete scheme and establish a result on the convergence of 
the discrete solution as the step size goes to zero. 

The Place of DMOC amongst Solution Methods for Optimal Control Problems. In Figure 3.1 
we present schematically different discretization strategies for optimal control problems: 

• In an indirect method, starting with an objective function and the Lagrange-d'Alembert principle we 
obtain via two variations (the first for the derivation of the Euler-Lagrange equations and the second for 
the derivation of the necessary optimality conditions) the Pontryagin maximum principle. The resulting 
boundary value problem is then solved numerically, e.g. by gradient methods (Cauchy [1847]; Kcllcy 
[I960]; Tollc [1975]; Bryson and Ho [1975]; Miclc [1980]; Chernousko and Luybushin [1982]), multiple 
shooting (Fox [I960]; Keller [1968]; Bulirsch [1971]; Deuflhard [1974]; Bock [1978]; Hiltmann [1990]) or 
collocation (Dickmanns and Well [1975]; Bar [1983]; Aschcr, Christiansen, and Russell [1979]). 

• In a direct approach, starting form the Euler-Lagrange equations we directly transform the problem into 
a restricted finite dimensional optimization problem by discretizing the differential equation. Common 
methods like e. g. shooting Kraft [1985], multiple shooting Bock and Plitt [1984], or collocation methods 
von Stryk [1993], rely on a direct integration of the associated ordinary differential equations or on its 
fulfillment at certain grid points (see also Betts [1998]; Pytlak [1999] for an overview of the current 
state of the art). The resulting finite dimensional nonlinear constrained optimization problem can 
be solved by standard nonlinear optimization techniques like sequential quadratic programming (Han 
[1976]; Powell [1978]). Implementations are found in software packages like DIRCOL von Stryk [2000], 
SOCS Betts and Huffmann [1998], or MUSCOD Lcinewcbcr [1999]). 

• In the DMOC approach, rather than discretizing the differential equations arising from the Lagrange- 
d'Alembert principle, we discretize in the earliest stage, namely already on the level of the variational 
principle. Then, we consider variations only on the discrete level to derive the restricted optimization 
problem and its necessary optimality conditions. 

This approach derived via the concept of discrete mechanics leads to a special discretization of the 
system equations based on variational integrators, which are dealt with in detail in Marsden and West 
[2001]. Thus, the discrete optimal control problem inherits special properties exhibited by variational 
integrators. In the following, we specify particular important properties and phenomena of variational 
integrators and try to translate their meaning into the optimal control context. 

Preservation of Momentum Maps. If the discrete system, obtained by applying variational integration 
to a mechanical system, inherits the same symmetry groups as the continuous system, the corresponding 
discrete momentum maps are preserved. For the forced case the same statement holds, if the forcing is 
orthogonal to the group action (see Theorem 2.2). 

On the one hand, this means for the optimal control problem, that if the control force is orthogonal to 
the group action, our discretization leads to a discrete system, for which the corresponding momentum map 
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Figure 3.1: Optimal control for mechanical systems: the order of variation and discretization for deriving 
the necessary optimality conditions. 



is preserved. On the other hand, in the case of the forcing not being orthogonal to the group action, the 
forced discrete Noether's theorem provides an exact coherence between the change in angular momentum 
and the applied control force via 

JV-l 

[ J C°( F id) ~ J L^}(lo,qi) ■£ = fd k (Qk,qk+i) ■ ^QxQ(qk,qk+i), 

k=0 

(see §4.2 for examples). 

Conservation of a Modified Energy. Variational integrators are symplectic, which implies that a certain 
modified energy is conserved (see for example Haircr, Lubich, and Wanner [2002]). This is an important 
property if the long time behavior of dynamical systems is considered. For the case of the optimal control 
of systems with long maneuver time such as low thrust space missions, it would therefore be interesting 
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to investigate the relation between a modified energy and the virtual work. However, this has not been 
considered within this paper. 



Implementation. Rather than using a configuration- momentum implementation of variational integrators 
as proposed in Marsden and West [2001], we stay on Q x Q. That means we just determine the optimal 
trajectory for the configuration and the control forces and reconstruct the corresponding momenta and 
velocities via the forced discrete Legendre transforms. This yields computational savings. A more detailed 
description of the computational savings compared to standard discretizations for optimal control problems 
is given in Remark 3.11. 



3.4 The Correspondence with Runge-Kutta Discretizations 

In this section we are going to show that the discretization derived via the discrete Lagrange-d'Alembert 
principle is equivalent to one resulting from a finite difference discretization of the associated Hamiltonian 
system via a symplectic partitioned Runge-Kutta scheme. 



Symplectic Partitioned Runge-Kutta Methods. As shown in Marsden and West [2001] (Theorem 
2.6.1), the discrete Hamiltonian map generated by the discrete Lagrangian is a symplectic partitioned Runge- 
Kutta method. As we will show, a similar statement is true for discrete Hamiltonian maps with forces. The 
resulting method is still a partitioned Runge-Kutta method, but no longer symplectic in the original sense 
since the symplectic form is not preserved anymore due to the presence of control forces. However, we still 
denote it as a symplectic method having in mind that the symplectic form is preserved only in absence of 
external control forces. 

A partitioned Runge-Kutta method for the regular forced Lagrangian system (L, is a map T*Q x 
jjs _^ t*q specified by coefficients bi, ay, bi, ay, i = l,...,s, and defined by (qo,Po,uo) i— > (qx, pi), where 



<7i = qo 



qo 



hy^bjQj, pi — po + h 2J bjPj , (3.25a) 

s s 

/i^JayQj, Pi = p Q + h2_jO-ijPj, (3.25b) 

j'=i i=i 
dL dL 

Pi = -griQuQi), Pi = -Q-(Qi,Qi) + h{Qi,Qi,Ui), (3.25c) 

i = 1, . . . , s, where the points (Qi, Pi) are known as the internal stages and C/j are the control samples given 
by Ui = uq{ = Ud(to + Cih). For ay = dy and bi = bi the partitioned Runge-Kutta method is a Runge-Kutta 
method. 

The method is symplectic (that is, it preserves the canonical symplectic form il on T*Q in absence of 
external forces) if the coefficients satisfy 

b.n.j ■ ';,</,, = bibj, i,j = 1, ...,s, (3.26a) 

bi = hi, i = l,...,s. (3.26b) 

Since discrete Lagrangian maps are symplectic, we can assume that we have coefficients satisfying (3.26) and 
write a discrete Lagrangian and discrete Lagrangian control forces that generate the corresponding symplectic 
partitioned Runge-Kutta method. Given points (go, Si) € Q x Q, we can regard (3.25) as implicitly defining 
POiPi, Qi, Pi, Qi and Pi for i = 1, . . . , s. Taking these to be defined as functions of (qo, q±), we construct a 
discrete Lagrangian 

s 

L d (q , qi ,h) = hJ2 b iL(Qi,Qi), (3.27) 
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and left and right discrete forces as 

S • dQi 

fd{qo,qi,u ,h) = hy bift{Qi,QuUi) ■ — jL , (3.28a) 

fd(lo,qi,u Q ,h) = h'S2b i f L {Q i ,Q i ,U i ) ■ (3.28b) 

For fixed uq the corresponding forced discrete Hamiltonian map is exactly the map (qo,po) > (qi,Pi) which 
is the symplectic partitioned Runge-Kutta method for the forced Hamiltonian system (2.6). 

Theorem 3.8. The discrete Hamiltonian map generated by the discrete Lagrangian (3.27) together with the 
discrete forces (3.28) is a partitioned Runge-Kutta method (which is symplectic in the unforced case). 

Proof. This is straight forward to the proof for unforced systems in Marsden and West [2001] (Theorem 
2.6.1). ■ 

Optimal Control and Runge-Kutta Discretizations. In this paragraph, we carry forward the results 
of the previous section into the context of optimal control problems. To formulate the optimal control problem 
for the discrete system in terms of Runge-Kutta discretizations, we have to give appropriate expressions for 
the boundary conditions, the path constraints, and the objective functional. Concerning the dynamics of 
the mechanical system, we have already seen that the discretization obtained via the discrete Lagrange- 
d'Alembert principle can be rewritten as a Runge-Kutta scheme for the corresponding mechanical system in 
terms of (p, q) as in (3.25). Since we consider regular Lagrangians and regular Hamiltonians, we reformulate 
(3.25c) with the help of the Hamiltonian as 

3U 

Qki -Qp{Qki,Pki), 

Pki = —-Q^{Qki,Pki) + fniQki, Pki,Uki), 

where the additional index k denotes the dependence of the intermediate state and control variables on the 
time interval [tfc,tfc+i]. 

Boundary Conditions. Due to a formulation of the optimal control problem within the Hamiltonian 
framework, we use the same formulation for the boundary constraint as the one for the continuous Hamilto- 
nian optimal control problem 3.3 evaluated at ((?at,pjv), which reads 

r(qN,PN,q T ,p T ) = 0. 

Discrete Path Constraints. Again we use the formulation of the Hamiltonian optimal control problem 
3.3 and enforce the path constraints to hold in each time step (Qi, Pi, C/j) as 

h(Qki,Pki,U ki )>0, k = 0,...,N-l, i = l,...,s. 

Discrete Objective Function. We construct the discrete objective function in the same way as the 
discrete Lagrangian (3.27). However, corresponding to the Hamiltonian formulation evaluated in each substep 
(Qi, Pi,Ui), it now reads 

JV-l s 

Jd(qd,Pd,Ud) = ^2 h^2biC(Qki,P ki ,U k i) + $(q N ,p N ), 

k=0 t=l 

where the final constraint holds for the node corresponding to the final time T. 
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Combining terms, the discrete optimal control problem from the Hamiltonian point of view reads 



min Jd(q<i,Pd,Ud) (3.29a) 
qd,Pd,fd 



subject to 



gfc+l = gfc + hj^bjQkj, Qki = qk + h^aijQkj, q = q°, (3.29b) 
3=1 3=1 

S S 

Pk+i=Pk + h^rbjP kj , P k j =Pk + h'Y]a i jPkj, Po = p°, (3.29c) 

3=1 3=1 

OH OH 

Q kl = -Q^{Qki,Pki), Pki = —g-(Qki,P k i) + fH(Qki,Pki,U ki ), (3.29d) 

Q< h(Qki,Pki,U ki ), k = Q,...,N-l,i = l,...,s, (3.29e) 
= r{q N ,p N ,q T ,p T ). (3.29f) 

Problem (3.29) is the finite dimensional optimization problem resulting from the Hamiltonian optimal control 
problem 3.3 that is discretized via a symplectic partitioned Runge-Kutta scheme. 

Theorem 3.9 (Equivalence). Given a discrete Lagrangian and discrete forces as defined in (3.27) and (3.28), 
respectively. Then, the discrete Lagrangian optimal control problem defined in (3.15) and the problem (3.29) 
resulting from discretizing the Hamiltonian optimal control problem by a symplectic partitioned Runge-Kutta 
scheme are equivalent in the sense that both problems have the same set of solutions. 

Proof. Assume (q d ,u d ) is a solution of (3.15). By using the discrete Legendre transform we obtain an 
optimal solution in terms of the discrete momenta as (q d ,p* d ,u* d ). To prove that {qd,Pdi u *d) * s a ^ so an 
optimal solution for problem (3.29) we have to show feasibility and optimality. Theorem 3.8 and appli- 
cation of the forced discrete Legendre transform to the boundary conditions give us the equivalence of 
equations (3.15b)-(3.15d),(3.15f) and (3.29b)-(3.29d), (3.29f). For the choice of the discretizations of the 
curves q(t) and u(t) the inequality condition (3.15e) reads as h{Q k i, Qk,i+i> U k i) > for i = 1, . . . , s — 1 and 
h(Qks,Qk+i.o,Uks) > 0. Again, due to the forced discrete Legendre transform h and h defined in (3.29c) 
determine the same solution set. This makes {q*dTP*dT u d) a feasible solution of the problem (3.29). Optimality 
holds due to Jd{q* d ,P*^ u d) = J d{q*d^* d ) ^ Jd{qd,u d ) = Jd(qd,Pd,u d ) for all feasible points (q d ,Pd) € T*C d (Q) 
and u d G Cd(U) (global optimality), or for all feasible points (qd,Pd, u d) in a neighborhood of {q* dl p d i u d) 
(local optimality). Analogous, we can show that an optimal solution (9j,p2' u d) of problem (3.29) is also an 
optimal solution for problem (3.15). ■ 



Mayer Formulation. Similar to what we did for the discrete Lagrangian optimal control problem, (3.29) 
can be transformed into an optimal control problem of Mayer type as follows: Analogous to what we did in 
equation (3.18), we introduce additional variables yd — {yi}fLo 85 

Vo = 0, 

l-l s 

V 1 = J2 h I2 bi C(Qki,Pki,Uki), l = l,...,N, 

k=0 i=l 

yielding the discrete optimal control problem of Mayer type as 

min $(qN,PN,yN) = Vn + ®(qN,PN) (3.30a) 

qd,Pd,u, d 

subject to (3.29b), (3.29c), (3.29d), (3.29e) and (3.29f). We obtain exactly the same problem by discretizing 
the continuous Hamiltonian optimal control problem of Mayer type with the same partitioned Runge-Kutta 
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discretization for the extended system of differential equations 

g(*) = *(«(*),#(*), u(i)), g(0)=g°, 
P(t)=r)(q(t),p(t),u(t)), p(0)=p°, 

with 9 - (<z,y), p = p, !/(«(*).*(*).«(*)) = (V p fl"(g(t),l>(t)),C'(«(t),p(t),u(t))) ) t/($(t),p(t),«(i)) - 
-Vpff( 9 (i),p(t))+/ H ( g (t) > p(t) J «(*)) 1 g° = (<?°,0) and p°=p°. 

Examples 3.10. 

(a) With 6=1 and a = | we obtain the implicit midpoint rule as a symplectic Runge-Kutta scheme, 
that is the partitioned scheme reduces to a standard one-stage Runge-Kutta scheme. The resulting 
discretization is equivalent to the discretization derived via the Lagrangian approach with midpoint 
quadrature. 

(b) The standard Lobatto IIIA-IIIB partitioned Runge-Kutta method is obtained by using the Lobatto 
quadrature for the discrete Lagrangian. 

Remark 3.11. 

1. The formulations of the discrete optimal control problem on the one hand based on the discrete 
Lagrange-d'Alembert principle and on the other hand based on a Runge-Kutta discretization of the 
Hamiltonian dynamics are equivalent in the sense that the same discrete solution set is described. 
Note however, that the Lagrangian formulation needs less discrete variables and less constraints for 
the optimization problem, as it is formulated on the configuration space only, rather than on the space 
consisting of configurations and momenta: For q £ W 1 and N intervals of discretization we obtain with 
the Lagrangian approach (As + l)n unknown configurations q%, k = 0, . . . , N — 1, v = 0, . . . , s with 
9a! = <7fc+n k — 1, . . . , A — 1 and nN(s — 1) + n(N — 1) extended discrete Euler-Lagrange equations, 
so altogether, (As — l)n constraints for the optimization problem excluding boundary conditions. The 
Runge-Kutta approach for the Hamiltonian system yields 2(As + l)n unknown configurations and 
momenta and 2nN + 2nN(s — 1) = 2nNs equality constraints. Thus, via the Runge-Kutta approach 
we obtain twice as many state variables and ( As + l)n more equality constraints such that the resulting 
optimization problem is numerically more expensive. Comparisons concerning the computational effort 
are presented in Section 4.2.2. 

2. Another advantage of the variational approach over the symplctic Runge-Kutta discretization is that 
one does not have to handle any conditions on the coefficients (such as condition (3.26)) to enforce 
symplecticity. The symplecticity of the DMOC discretization results naturally from the variational 
structure of this approach. 

3.5 The Adjoint Systems 

The adjoint system provides necessary optimality conditions for a given optimal control problem. In the 
continuous case the adjoint system is derived via the Pontryagin maximum principle. The KKT equations 
provide the adjoint system for the discrete optimal control problem. Hagcr [2000] derives a transformed 
adjoint system using standard Runge-Kutta discretizations. He identifies order conditions on the coefficients 
of the adjoint scheme up to order 4. Bonnans and Laurent-Varin [2006] extend these up to order 7. 

By using the same strategy as Hager in Hagcr [2000] , we here show that DMOC leads to a discretization 
of the same order for the adjoint system as for the state system. Therefore no additional order conditions on 
the coefficients are necessary. In the following we ignore path and control constraints and restrict ourselves 
to the case of unconstrained optimal control problems. 

Continuous Setting. We consider a Hamiltonian optimal control problem in Mayer form without final 
constraint, path and control constraints: 
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Problem 3.12. 



subject to 



min$(g(T),p(T)) (3.31a) 



q.p.u. 



q = v{q,P), (3.31b) 
P=v(l,P,u), (3.31c) 
(g(0),p(0)) = (q°,p°), (3.31d) 



with q,p e W^ftO,!],®"), U S L°°([0,r|,R m ), 1/(9, p) - V p H(q,p), V (q,p,u) = -V q H(q,p) + f B (q,p, 



ii 



We also write x(t) = (q(t),p(t)), x° = (q°,p°) and 

f(x, u) = 



v{q,p) 
rj(q,p,u) 

such that (3.31b)- (3.31d) reads as x — f(x,u),x(0) = x°. 

Along the lines of Hagcr [2000], we now formulate the assumptions that are employed in the analysis 
of DMOC discretizations of Problem 3.12. First, a smoothness assumptions is required to ensure regularity 
of the solution and the problem functions. Second, we enforce a growth condition that allows for having a 
unique solution for the control function of the optimal control problem. 

Assumption 3.13 (Smoothness). For some integer k>2, Problem 3.12 has a local solution (x*,u*) which 
lies in W K '°° x W^' 1 ' 00 . There exists an open set fl C M. 2n x W n and p > such that 2 B p (x*(t),u*(t)) C il 
for every t G [0, T] . The first k derivatives of v and r] are Lipschitz continuous in £1, and the first k derivatives 
of & are Lipschitz continuous in B p {x*{T)). 

First order optimality conditions. Under Assumption 3.13 there exits an associated Lagrange multiplier 
ip* = (t/; 9 '* , ipP'*) £ W K '°° for which the following form of the first-order optimality conditions derived via 
the Pontryagin maximum principle is satisfied at (x* , ip* , u*): 

q=u(q,p), g(0) =9°, (3.32a) 

p = T!(q,p,u), p(0)=p°, (3.32b) 

r = -V q H(q,p, u), r(T) = V 9 <%(T),p(T)) (3.32c) 

r = ~V p H(q lP ^\r,u), V P (T) = V p $(q(T),p(T)) (3.32d) 

V tt W(g(t), p(t), il> q (t), V p (t), «(t)) = for alH e [0,T], (3.32e) 

with the Hamiltonian TL defined by 

H(q,p, r,V,u) = ^v{q,p)+rv{q,P,u), (3.33) 

where ip q and ip p are row vectors in R". We also write H.(x,ip,u) := ipf(x,u) such that (3.32c) and (3.32d) 
read as 

V> = -V x H(x, i>, u), i>{T) = V x $(x(T)). 

Second order optimality conditions. In order to formulate the second-order sufficient optimality con- 
ditions we define the following matrices (cf. Hagcr [2000]): 



A(t) 


= V x f(x*(t) >U *(t)), 




B(t) 


= Vj(x*(t), u*(t)), 




V(t) 


= V 2 4>(x*(T)), 




P(t) 


= V m H(x*(t) ) il> m (t), 


«*(*)), 


R(t) 


= V M «(i*(t),f(t), 


u*(t)), 


S(t) 


= V xu H(x*(t),t/i*{t), 


«*(*))• 



2 Bp(z) is the closed ball centered at z with radius p. 
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Let B be the quadratic form defined by 




where (•, •) denotes the usual L 2 



inner product. 



Assumption 3.14 (Coercivity). There exists a constant a > such that 



B(x,u) > allulll, for all (x,u) £ M.. 



where M = {(x,u) £ H 1 x L 2 \ x(0) = and x = Ax + Bu for a.e. t £ [0,T]}. 

Under the assumptions Smoothness and Coercivity, the control uniqueness property holds (cf. Hager 
[1979]): the Hamiltonian TL has a locally unique minimizcr u* = u(x, tp) in the control, depending Lipschitz 
continuously on x and tp. Let (j> = (4> q > 4^) denote the function defined by 



Additionally, let rj(x,tp) = tli. x i u *)- By substituting u* = u(x,iJj) into (3.32) one obtains a two-point 
boundary-value problem 



Discrete Setting. It was shown in §3.4 that the discrete Lagrangian control problem is equivalent to 
Problem (3.29) (resp. (3.30)). Using some transformation and change of variables (see Hager [2000] and 
Ober-Blobaum [2008]), by reversing the order of time and applying the control uniqueness property, we 
obtain the following version of the first-order necessary optimality conditions associated with (3.29) (the 
Karush-Kuhn- Tucker conditions) : 



4>"{x^) 



X7 q H(x,ip,u*) 7 
W p H(x, -0, u*). 



p = v(q,p, 4> q ,V), 
y> = ^(<7,p,^,V p ), 



p(Q)=p°, 
0«(T) = V,$(g(T),p(T)), 
r(T)=V p $(q(T),p(T)). 



(3.34a) 
(3.34b) 
(3.34c) 
(3.34d) 
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Qk+i = Qk + h^2 / b i v{Q ki , P ki ), q = q°, (3.35a) 



i=l 



Pfc+i=p fc + /iy]MOw.flw»*L*w)» Po=P°, (3.35b) 



i=l 



Qfei = qk + h^a q i:j iy(Qkj,Pkj), (3.35c) 



S 

P kl =Pk + hY / a p ijV (Q kj , P kj , n 3 , *lj ) > (3-35d) 

3=1 
s 

^Ux = i>i + hY,h4> q {Qki,Pki^li^lil 4>% = V 9 *(qN,PN), (3.35e) 
1=1 



^ + i = V£ + &^M p (Qki,fiw, ^ = V P $(^, PJV ), (3.35f) 

i=l 
s 

Hi = K + htty 9 (QkS,Pk i ,* q kj ,* p kj )> (3-35g) 



* ? « = ^ + '«L V(Q fci , i^, *U ( 3 - 35h ) 

3=1 

a|,- = . J Jt , a?.= , J • (3.35i) 



Assuming h small enough and x k — (qk>Pk) near x*{t k ) = (q* (tk) , P* (tk)) the intermediate states Q k i and 
Pfci and costates ^ q ki and are uniquely determined since u(x, tp) depends Lipschitz continuously on x near 
x*(t) and ip near ip*(t) for any £ € [0,T]. This follows from smoothness and the implicit function theorem 
as for standard Runge-Kutta discretization as stated in Hagcr [2000] (sec for example Flatto [1976]). Thus, 
there exists a locally unique solution (Q k i, P k i, ^kv^kd' 1 < * < s °f (3.35e)-(3.35h). 

The scheme (3.35) can be viewed as a discretization of the two-point boundary- value problem (3.34). 
To ensure a desired order of approximation for this two-point boundary-value problem, Hagcr [2000] derives 
order conditions on the coefficients of the Runge-Kutta scheme via Taylor expansions. In our case, however, 
we more easily obtain the following result concerning the order of consistency. 

Theorem 3.15 (Order of consistency). If the symplectic partitioned Runge-Kutta discretization of the state 
system is of order k and bi > for each i, then the scheme for the adjoint system is again a symplectic 
partitioned Runge-Kutta scheme of the same order (in particular we obtain the same schemes for (q,p) and 

Proof. Starting with a symplectic partitioned Runge-Kutta discretization for the state system (3.31), the 
discrete necessary optimality conditions are given as the adjoint system (3.35e-3.35i). This system is again 
a partitioned Runge-Kutta scheme for the adjoint equations with coefficents &j,af-,a?- . Substituting the 

symplecticity condition aP- — b ' bj ^ : ' ari into equation (3.35i), the coefficients are determined as 

-q p 

a ij 

-P Q 

<i = 4- 

Since the coefficients of the Runge-Kutta scheme of the adjoint system for (ip p ,ip q ) are the same as the 
coefficients of the Runge-Kutta scheme of the state system for (q,p) defined in (3.35a-3.35d), the adjoint 
scheme is of same order. ■ 

With Theorem 2.5 and 3.8 one obtains the following: 
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Lemma 3.16. Let a Lagrangian optimal control problem with regular Lagrangian L and Lagrangian control 
force Jl be given. If the discrete Lagrangian Ld (3.27) and the discrete control force f d (3.28) with bi > 
for each i are both of order k, then the corresponding adjoint scheme is also of order k. 

3.6 Convergence 

The purpose of this section is to establish the convergence of solutions of the discrete optimal control problem 
(3.15) to a solution of the Lagrangian optimal control problem (3.3) as the step size h goes to 0. As in Section 
3.5, we restrict ourselves to the case without path and control constraints and without final point constraint. 
Our convergence statement is a direct application of the one in Dontchev, Hager, and Veliov [2000] (cf. 
also Hager [2000] and Obcr-Blobaum [2008]). In the previous section, we derived the adjoint system and 
determined the order of consistency. A standard strategy to show convergence of a discrete scheme is to 
prove consistency and stability. This is abstractly stated within the following result: 

Theorem 3.17. Hager [2001] Let X be a Banach space and let y be a linear normed space with the norm 
in both spaces denoted by \\ • ||. Let T : X i— > 2^ be a set-valued map, let C : X i— ► y be a bounded, linear 
operator, and let T : X i— > y with T continuously Frech'et differentiate in B r (w*) for some w* G X and 
r > 0. Suppose that the following conditions hold for some 5 G y and scalars e, A and a > 0: 

(PI) T(w*) + 5e T{w*). 

(P2) \\VT(w) - C\\ < e for all w G B r (w*). 

(P3) The map (J-—C)^ 1 is single-valued and Lipschitz continuous in B a (ir),ir — (T —jC)(w*), with Lipschitz 
constant A. 

IfeX < I, er < a, and \\6\\ < (1 
Moreover, we have the estimate 



Consistency corresponds to assumption (PI) and the bounds on the norm of S, stability corresponds to 
assumption (P3) and the bound on the Lipschitz constant A for the linearization, and convergence is stated 
in (3.36). 

The following convergence result is formulated in terms of the averaged modulus of smoothness of the 
optimal control. If J C M is an interval and v : J — > M n , let co(v, J; t, h) denote the modulus of continuity: 

u>(v,J;t,h)=swp{\v(s 1 )-v(s 2 )\ : si, s 2 G [t - h/2, t + ft/2] n J}. (3.37) 
The averaged modulus of smoothness r of v over [0, T] is the integral of the modulus of continuity: 

T (v;h) = [ w(u,[0,T];£,ft)dt 
Jo 

It is shown in Scndov and Popov [1988] that lim^o t(v; ft) = if and only if v is Ricmann integrable, and 
t(v; ft) < eft with constant c if v has bounded variation. 

Theorem 3.18 (Convergence). Let (x*,u*) be a solution of the Lagrangian optimal control problem (3.3a)- 
(3.3d). If smoothness and coercivity hold, the discrete Lagrangian Ld (3.27) and the discrete control force 
fj 1 (3.28) are of order k with bi > for each i, then for all sufficiently small ft there exists a strict 
local minimizer (x h ,u h ) of the discrete Lagrangian optimal control problem (3.15) and an associated adjoint 
variable ip satisfying the first-order necessary optimality conditions such that 

m^Jx h k -x*(t k )\ + \^-r(tk)\ + \u(xl^)-u*(t k )\ Kch*- 1 U + rt^u*;^^ , (3.38) 

where u(x k ,ip k ) is a local minimizer of the Hamiltonian (3.33) corresponding to x — x% and ifj = i/) k . 



\e)r/\, then there exists a unique w E B r (w*) such that T(w) E J-(w). 
||«,_t U *||< T A-|| < y||. (3.36) 
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Proof. By Theorems 2.5 and 3.8 we know that a discrete Langrangian and discrete forces both of order k lead 
to a symplectic partitioned Runge-Kutta discretization of order k for the state system. Because of smoothness 
and coercivity we can build up a discrete adjoint scheme with eliminated control that approximates the 
continuous adjoint scheme with order k (see Lemma 3.16). This leads in Hager's terminology to a Runge- 
Kutta scheme of order k for optimal control, and therefore Hager's convergence result for standard Runge- 
Kutta schemes in Hagcr [2000], Theorem 2.1, is directly applicable. ■ 

Remark. Note that the estimate for the error in the discrete control in (3.38) is expressed in terms of 
u(x^,tp^) not Uk.. This is due to the fact, that we derive the estimate via the transformed adjoint system 
with removed control due to the control uniqueness property. In Dontchcv, Hager, and Vcliov [2000] the 
estimate is proved in terms of for Runge-Kutta discretization of second order. 

Remark. Theorem 3.18 can be extended to optimal control problems with constraints on the control 
function u(t) as it was done in Hagcr [2000] for Runge-Kutta discretizations of order 2. 

4 Implementation and Applications 
4.1 Implementation 

As a balance between accuracy and efficiency we employ the midpoint rule for approximating the relevant 
integrals for the example computations in the following section, that is we set 



k = 0, . . . ,N - 1. Here, = f£ = \ f L [ qk+x 2 +qk , qk+1 h ^ ,u k+1/2 ) are used as the left and right discrete 
forces with q). — q{tk) and u k+ i/ 2 = U ( tfc+ * fc+1 J . 



SQP Method. We solve the resulting finite dimensional constrained optimization problem by a standard 
SQP method as implemented for example in the routine fmincon of MATLAB. For more complex problems 
we use the routine nag_opt jnlp_sparse of the NAG library 3 . Since SQP is a local method, different initial 
guesses can lead to different solutions. This has been observed in almost all our example computations. 

Automatic Differentiation. To compute an optimal solution of the optimization problem the SQP 
method makes use of the first and second derivatives of the constraints and the objective function. In 
the case where no derivatives are provided, the algorithms that are used approximate those by finite dif- 
ferences. This approximation is time-consuming and round-off errors in the discretization process occur 
leading to worse convergence behavior of the algorithm. To avoid these drawbacks we make use of the 
concept of Automatic Differentiation (AD) (see Gricwank [2000]; Rail [1981]; Wengcrt [1964]), a method to 
numerically evaluate the derivative of a function specified by a computer program. AD exploits the fact 
that any computer program that implements a vector function y = F(x) (generally) can be decomposed into 
a sequence of elementary assignments, any one of which may be trivially differentiated by a simple table 
lookup. These elemental partial derivatives, evaluated at a particular argument, are combined in accordance 
with the chain rule from differential calculus to yield information on the derivative of F (such as gradients, 




www. nag . com 
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tangents, and the Jacobian matrix). This process yields (to numerical accuracy) exact derivatives. Because 
symbolic transformations occur only at the most basic level, AD avoids the computational problems inherent 
to complex symbolic computation. In particular, for some of our optimal control problems, we implemented 
the package ADOL-C (Automatic Differentiation by OverLoading in C++, Walthcr, Kowarz, and Gricwank 
[1996]) that has been written primarily for the evaluation of gradient vectors (rows or columns of Jacobians). 
It turns out that the optimization process performs in a faster and more robust way when providing the 
derivatives by automatic differentiation rather than via finite differences. 



4.2 Applications 

In this section, we numerically compare DMOC to a collocation method by means of two problems: a low 
thrust orbital transfer and the optimal control of a two-link manipulator. For this comparison, we apply a 
collocation method of order 2 to two different models: the Hamiltonian system with coordinates (q,p) as 
well as the system formulated on tangent space with coordinates (q,q) — (q,v). 



4.2.1 Low Thrust Orbital Transfer 



In this application we investigate the problem of optimally transferring a satellite with a continuously acting 
propulsion system from one circular orbit around the Earth to another one. 



Model. Consider a satellite with mass m which moves in the gravitational field of the Earth (mass M). 
The satellite is to be transferred from one circular orbit to one in the same plane with a larger radius, while 
the number of revolutions around the Earth during the transfer process is fixed. In 2d-polar coordinates 
q = (r, ip), the Lagrangian of the system has the form 

t/ -\ 1 i-i 2 -2\ Mm 

HQ, l) = o m ( r + r <P )+ 7 > 

2 r 

where 7 denotes the gravitational constant. Assume that the propulsion system continuously exhibits a force 
u in the direction of motion of the satellite, such that the corresponding Lagrangian control force is given by 




Boundary Conditions. Assume further that the satellite initially moves on a circular orbit of radius 
r°. Let (r(0),yj(0)) = (r°,0) be its position at t — 0, then its initial velocity is given by f(0) = and 
<p(0) — y/^/M/lr ) 3 . Using its thruster, the satellite is required to reach the point (r T ,0) at time T = 

8-yM ( r ° rT ) 3 an d, without any further control input, to continue to move on the circle with radius r T . 
Here, d is a prescribed number of revolutions around the Earth. Thus, the boundary values at t — T are 
given by (r(T),tp(T)) = (r T , 0) and (f(T), ip(T)) = (0, ^/ 1 M/{r T f). 



Objective Functional. During this transfer, our goal is to minimize the control effort, correspondingly 
the objective functional is given by 



J(<7, q, u) = / u(t) dt 



T 

2 



Results. We compute the transfer from an orbit of radius r° = 30 km to one of radius r T = 330 km 
around the Earth. First, we investigate how well the balance between the change in angular momentum and 
the amount of the control force is preserved. Due to the invariance of the Lagrangian under the rotation if, 
the angular momentum of the satellite is preserved in the absence of external forces (as stated in Noether's 
theorem). However, in the presence of control forces, equation (2.7) gives a relation between the forces and 
the evolution of the angular momentum. 
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log 2 (N=2 k ) 

Figure 4.1: Comparison of the accuracy of the computed open loop control for DMOC and a collocation 
approach: difference of force and change in angular momentum in dependence on the number of discretization 
points. 



In Figure 4.1, we compare the amount of the acting force with the change in angular momentum 
in each time interval. For the solution resulting from DMOC and the collocation approach applied to 
the Hamiltonian system, the change in angular momentum exactly equals the sum of the applied control 
forces (to numerical accuracy) . The collocation method of second order corresponds to a discretization via 
the implicit midpoint rule. Thus, the optimization problem resulting from DMOC is equivalent to that 
obtained by applying collocation to the Hamiltonian formulation of the system as shown in Theorem 3.9 
and Example 3.10a). Obviously, we obtain equal solutions by applying both methods. These results are 
consistent with the well-known conservation properties of variational integrators, that provide discretizations 
that preserve the continuous properties as momentum maps in the discrete setting in a natural way. On 
the other hand, the collocation method applied to the tangent space system described in velocities fails to 
capture the change in angular momentum accurately because the discrete tangent space formulation destroys 
the discrete Hamiltonian structure and the resulting (also unforced) scheme is not momentum-preserving 
anymore. 

In Figure 4.2 we show the convergence rates for all three methods. A reference trajectory is computed 
with N = 1024 discretizations points and time step h — 2.9- 10 -3 . The error in the configuration and control 
parameter of the discrete solution with respect to the reference solution is computed as maxto,...,jv Is(ife) — 
<7ref(*fc)| an d niaxfc = o ) 7v \u(tk)— u [c { (ife)|, respectively, where |-| is the Euclidean norm. For all three methods 
the convergence rate for the configuration and control trajectory is 0(h 2 ), as predicted by Theorem 3.18. 

Still, DMOC is advantageous regarding the computational efficiency. Due to its formulation on the 
discrete configuration space, DMOC only uses I w 0.6 of the number of variables of the collocation approach 
(cf. Remark 3.11). Figure 4.3 shows the number of all iterations that the SQP solver performs in order to 
solve the quadratic subproblems (minor iterations). We observe that for DMOC the SQP solver needs about 
1.7 times less iterations in comparison to collocation. This is in exact aggreement with the reduced number 
of variables. 

4.2.2 Two-link Manipulator 

As a second numerical example we consider the optimal control of a two-link manipulator. Again, we compare 
DMOC to a collocation method of the same order. 



4.2 Applications 33 




'l 2 3 4 5 6 2 1 2 3 4 5 6 

-log(h) -iog(h) 



Figure 4.2: Approximation error for the states (left) and the controls (right) in dependence of the step size 
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Figure 4.3: Comparison of the number of iteration steps performed by the SQP solver for DMOC and a 
collocation approach in dependence on the number of discretization points. 
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Model. The two-link manipulator (see Figure 4.4) consists of two coupled planar rigid bodies with mass 
rrii, length lj and moment of inertia Ji, i — 1,2, respectively. For j £ 1,2, we let denote the orientation of 
the «th link measured counterclockwise from the positive horizontal axis. If we assume one end of the first 
link to be fixed in an inertial reference frame, the configuration of the system is specified by q = {0\, 9 2 ). 




2 



Tl 



Figure 4.4: Model of a two-link manipulator. 

The Lagrangian is given via the kinetic and potential energy 

K(q, q) = i(mi + Am 2 )l\el + ^m 2 l 2 6 2 + l -m 2 hh cos {9 1 - d 2 )B x B 2 + X -J x d\ + l -J 2 B 2 

and 

V(q) = 2 TO i9'i sill ^i + m 29k sin 6 X + ^m 2 gl 2 8 2 , 

with the gravitational acceleration g. Control torques ri, t 2 are applied at the base of the first link and at 
the joint between the two links. This leads to the Lagrangian control force 

to.*)- (V)- 

Boundary Conditions. The two-link manipulator is to be steered from the stable equilibrium point 
q° = (—§,—§) with zero angular velocity q° — (0,0) to the unstable equilibrium point q T = (~, |) with 
velocity q T = (0,0). 

Objective Functional. For the motion of the manipulator the control effort 

J(n,r 2 )= jT T i(7?(t) + 7f(t))dt 

is to be minimized, where we fix the final time T — 1. 

Results. In Figure 4.5 we show a) the resulting cost and b) the difference of the amount of force (including 
the control and the gravitational force) and the change in angular momentum in dependence on the number 
of discretization points for all three methods. As expected, we obtain (numerically) identical solutions for 
DMOC and the equivalent collocation method for the Hamiltonian formulation. The midpoint rule applied to 
the tangent space formulation performs equally well with respect to the objective value evolution. However, 
as in the previous example it does not reflect the momentum-force consistency as good as the other methods 
as shown in Figure 4.5. In Figure 4.6 the convergence rates are depicted. Here, a reference trajectory has 
been computed with N — 512 discretizations points and time step h = 1.9 • 10 -3 . For all three methods the 
convergence rate for the configuration and control trajectory is 0(h 2 ), as expected for a scheme of second 
order accuracy. 
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a) b) 

Figure 4.5: Comparison of the accuracy of the computed open loop control for DMOC and a collocation 
approach: a) approximated cost b) difference of force and change in angular momentum in dependence on 
the number of discretization points. 



In Figure 4.3 the number of minor iterations that the SQP solver performs to solve the quadratic 
subproblems is shown. Similar to what we have seen in the orbital transfer, we observe that the SQP solver 
needs about 1.5 times less iterations when using DMOC rather than a collocation approach. This factor 
reflects the fact, that DMOC uses only | ss 0.67 of the number of variables in comparision to the collocation 
approach. 

5 Conclusions and Future Directions 

In this paper we developed a fully discrete formulation of the optimal control of mechanical systems. Based on 
discrete variational principles, the discrete optimal control problem yields a strategy, denoted by DMOC, for 
the efficient solution of these kinds of problems. The benefits of using DMOC are twofold: On the one hand, 
in the presence of a symmetry group in the continuous dynamical system, it is known that the momentum 
map changes in a definite way. The use of discrete mechanics allows one to find an exact counterpart to 
this on the discrete level. In this paper, this behavior was shown numerically for specific examples. On 
the other hand, due to the fact that DMOC is implemented on the configuration level rather than on 
the configuration-momentum or configuration-velocity level, one gets significant computational savings for 
the number of iterations steps the SQP solver needs to solve the optimization problem. Again, this was 
demonstrated numerically in examples. 

Although the developed method works very successfully in many examples considered in this paper and 
others (see Junge, Marsdcn, and Obcr-Blobaum [2006]; Jungc and Obcr-Blobaum [2005]; Junge, Marsdcn, 
and Obcr-Blobaum [2005]; Lcycndcckcr, Obcr-Blobaum, and Marsdcn [2008]; Obcr-Blobaum [2008]), there 
are still some open challenges that will be investigated in future work. 

Global Minima One challenge is the issue of local versus global minima. The use of an SQP solver for 
solving the resulting optimization problem restricts one to local minima that are dependent on the specific 
initial guess used. One approach to overcome this issue is the use of DMOC primitives (Kobilarov [2008]; 
Frazzoli, Dahlch, and Feron [2005]). Here, the authors create a roadmap of feasible trajectories given by 
a graph where the edges represent small pre-computed DMOC segments referred to as DMOC primitives. 
A feasible solution of the optimal control problem corresponds to a specific concatenation of the DMOC 
primitives respecting the dynamics. The global optimal control can now be approximated by determining 
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Figure 4.7: Comparison of the number of iteration steps performed by the SQP solver for DMOC and a 
collocation approach in dependence on the number of discretization points. 
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the optimal path in the graph with a given cost function with a global search using techniques from dynamic 
programming. Initial demonstrations of this idea show that it can be carried out in real time. This, together 
with the fact that real dynamics is included via DMOC primitives are the significant advantages of this 
approach. 

Adaptive time-stepping For many applications, an adaptive time-stepping strategy is essential. For 
example, for problems in space mission design, the planned trajectories require a finer time-stepping nearby 
planets due to the strong influence of gravity, while for a transfer in nearly free space only few discretization 
points are necessary to accurately reflect the dynamics of the system. Here, different strategies such as error 
control based on the discretization grid under consideration and variational approaches could be investigated. 
For the variational approach a constraint is included to the Lagrange-d'Alembert principle that ensures time 
step control directly at the level of the discrete action (Kharcvych, Mullen, Lcyendeckcr, Tong, Marsdcn, 
and Desbrun [2008]). According to different adaption schemes different constraints can be included such as 
adaption to acceleration or the strenght of control forces. 

Miscellaneous The framework could be extended to the optimal control of mechanical systems with 
stochastic influence or contact problems, respectively making use of the theory of stochastic (Bou-Rabcc 
and Owhadi [2007]) and nonsmooth variational integrators (Fetecau, Marsden, Ortiz, and West [2003]), 
respectively. Due to the variational formulation of DMOC, an extension towards the optimal control of 
partial differential equations for the treatment of fluid and continuum mechanics might be interesting as well 
(see Marsden and Shkollcr [1999]; Marsden, Patrick, and Shkoller [1998] for basic extensions of variational 
integrators to the context of PDEs) . 
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